Evapotranspiration.f90 Source File
Compute evapotranspiration
Source Code
!! Compute evapotranspiration !|author: <a href="mailto:giovanni.ravazzani@polimi.it">Giovanni Ravazzani</a> ! license: <a href="http://www.gnu.org/licenses/">GPL</a> ! !### History ! ! current version 1.2 - 8th March 2023 ! ! | version | date | comment | ! |----------|-------------|----------| ! | 1.0 | 08/Jun/2021 | Original code | ! | 1.1 | 27/Jan/2023 | subroutine HargreavesSamaniModified added | ! | 1.2 | 08/Mar/2023 | etEnergyBalanceGround | ! !### License ! license: GNU GPL <http://www.gnu.org/licenses/> ! !### Module Description ! Module to compute evapotranspiration ! MODULE Evapotranspiration ! Modules used: USE DataTypeSizes, ONLY : & ! Imported Parameters: float, short USE Loglib, ONLY : & !Imported routines: Catch USE StringManipulation, ONLY : & !imported routines ToString USE GridLib, ONLY : & !Imported types: grid_real, grid_integer, & !Imported routines: NewGrid , GetDtGrid !, & !GridDestroy, ExportGrid, & !SetLongName, SetStandardName, & !SetUnits, SetCurrentTime, & !Imported parameters: !NET_CDF, ESRI_ASCII, & !ESRI_BINARY USE GridOperations, ONLY : & !imported routines: CRSisEqual, GridByIni, & !imported operators: ASSIGNMENT (=) USE IniLib, ONLY : & !Imported types: IniList, & !Imported routines: IniOpen, IniClose, & IniReadInt, IniReadReal, & KeyIsPresent, IniReadString, & SectionIsPresent USE Chronos, ONLY : & !imported definitions: DateTime, & !Imported routines: DayOfYear, & !Imported operations: OPERATOR (+), & OPERATOR (==), & OPERATOR (<), & ASSIGNMENT (=) USE Units, ONLY : & !imported parameters: pi, millimeter, day USE DomainProperties, ONLY : & !imported variables: mask, latCentroid USE MorphologicalProperties, ONLY : & !imported variables: dem, dem_loaded USE AirTemperature, ONLY : & !Imported variables: dtTemperature, temperatureAir USE AirTemperatureDailyMax , ONLY : & !imported variables: dtTemperatureDailyMax, temperatureAirDailyMax USE AirTemperatureDailyMin , ONLY : & !imported variables: dtTemperatureDailyMin, temperatureAirDailyMin USE AirTemperatureDailyMean , ONLY : & !imported variables: dtTemperatureDailyMean, temperatureAirDailyMean USE AirRelativeHumidity, ONLY : & !Imported variables: dtRelHumidity, relHumidityAir, hygrometers USE WindFlux, ONLY : & !Imported variables: dtWindSpeed, windSpeed, & anemometersWS USE SolarRadiation, ONLY : & !Imported routines: dtRadiation, netRadiation, netRadiationFAO USE Plants, ONLY : & !Imported variables: fvcover, fvcoverLoaded, & lai, laiLoaded, & plantsHeight, plantsHeightLoaded, & rsMinLoaded, rsMin USE ObservationalNetworks, ONLY : & !Imported types: ObservationalNetwork, & !Imported routines: ReadMetadata, ReadData !, & !CopyNetwork , & !Imported operators: !ASSIGNMENT (=) USE Utilities, ONLY : & !imported routines: GetUnit USE MeteoUtilities, ONLY: & !Imported routines: ReadField IMPLICIT NONE ! Global (i.e. public) Declarations: TYPE (grid_real) :: pet !!potential evapotranspiration rate [m/s] TYPE (grid_real) :: pe !!potential evaporation rate [m/s] TYPE (grid_real) :: pt !!potential transpiration rate [m/s] INTEGER (KIND = short) :: dtET !! time step computation [s] !Private declarations INTEGER (KIND = short), PRIVATE, PARAMETER :: ET_MODELS = 6 INTEGER (KIND = short), PRIVATE, PARAMETER :: PENMAN_MONTEITH = 1 INTEGER (KIND = short), PRIVATE, PARAMETER :: PRIESTLEY_TAYLOR = 2 INTEGER (KIND = short), PRIVATE, PARAMETER :: HARGREAVES = 3 INTEGER (KIND = short), PRIVATE, PARAMETER :: HARGREAVES_MOD = 4 INTEGER (KIND = short), PRIVATE, PARAMETER :: FAO56_PENMAN_MONTEITH = 5 INTEGER (KIND = short), PRIVATE, PARAMETER :: ENERGY_BALANCE = 6 TYPE (grid_integer), PRIVATE :: model_map TYPE (grid_integer), PRIVATE :: kc_code_map !!map of crop coefficient code TYPE (grid_real), PRIVATE :: kc_map !! crop coefficient map INTEGER (KIND = short), PRIVATE :: model_vector ( ET_MODELS ) !!defines active models INTEGER (KIND = short), PRIVATE :: model !!model to compute ET INTEGER (KIND = short), PRIVATE :: model_assignment !!method to assign computation !!1 = one method for the entire domain, 2 = a map with model codes INTEGER (KIND = short), PRIVATE :: fileunitKc INTEGER (KIND = short), PRIVATE :: dtKc = 0 INTEGER (KIND = short), PRIVATE :: interpolationMethodKc INTEGER (KIND = short), PRIVATE :: dtGridKc CHARACTER (LEN = 300), PRIVATE :: filenameKc TYPE (DateTime), PRIVATE :: tnewET TYPE (DateTime), PRIVATE :: tNewKc LOGICAL, PRIVATE :: compute_hargreaves = .FALSE. LOGICAL, PRIVATE :: useCropCoefficient = .FALSE. TYPE (ObservationalNetwork), PRIVATE :: kc_stations !! pseudo station to read kc time series values !Public routines PUBLIC :: ETinit PUBLIC :: PETupdate PUBLIC :: etEnergyBalanceGround !Local routines PRIVATE :: HargreavesSamani PRIVATE :: HargreavesSamaniModified PRIVATE :: PriestleyTaylor PRIVATE :: PenmanMonteith PRIVATE :: FAO56PenmanMonteith PRIVATE :: KcUpdate !======= CONTAINS !======= ! Define procedures contained in this module. !============================================================================== !| Description: ! Initialize evapotranspiration computation SUBROUTINE ETinit & ! ( inifile, time ) IMPLICIT NONE !Arguments with intent(in): CHARACTER (LEN = *), INTENT(IN) :: inifile TYPE (DateTime), INTENT(IN) :: time !local declarations: TYPE (IniList) :: etini INTEGER (KIND = short) :: i, j !------------------------------end of declarations----------------------------- !open and read configuration file CALL IniOpen (inifile, etini) !dt IF (KeyIsPresent('dt', etini) ) THEN dtET = IniReadInt ('dt', etini) ELSE CALL Catch ('error', 'Evapotranspiration', & 'dt missing in configuration file' ) END IF !initialize PET map CALL NewGrid (pet, mask) CALL NewGrid (pe, mask) CALL NewGrid (pt, mask) !model assignment method IF (KeyIsPresent('model-assignment', etini ) ) THEN model_assignment = IniReadInt ('model-assignment', etini) ELSE CALL Catch ('error', 'Evapotranspiration', & 'model-assignment missing in configuration file') END IF !set model IF (model_assignment == 1) THEN model = IniReadInt ('model', etini ) CALL NewGrid (model_map, mask) CALL CheckSpecificProperties (model, etini) model_map = model ELSE model = -1 !check map is present IF ( .NOT. SectionIsPresent ('model-map', etini) ) THEN CALL Catch ('error', 'Evapotranspiration', & 'model-map missing in configuration file') END IF !read map CALL GridByIni (etini, model_map, 'model-map') IF ( .NOT. CRSisEqual (mask, model_map, .TRUE.) ) THEN CALL Catch ('error', 'Evapotranspiration', & 'something seems wrong in model-map') END IF !scan for et models model_vector = -1 DO i = 1, model_map % idim DO j = 1, model_map % jdim IF (model_map % mat(i,j) /= model_map % nodata) THEN model_vector (model_map % mat (i,j)) = model_map % mat (i,j) END IF END DO END DO DO i = 1, ET_MODELS CALL CheckSpecificProperties (model_vector (i), etini) END DO END IF !crop coefficient IF ( KeyIsPresent ('use-crop-coefficient', etini) ) THEN IF (IniReadInt ('use-crop-coefficient', etini) == 1 ) THEN useCropCoefficient = .TRUE. END IF END IF IF (useCropCoefficient) THEN IF ( .NOT. SectionIsPresent ('crop-coefficient', etini) ) THEN CALL Catch ('error', 'Evapotranspiration', & 'crop-coefficient section is missing in configuration file' ) END IF !allocate map CALL NewGrid (kc_map, mask) !file filenameKc = IniReadString ('file', etini, section = 'crop-coefficient') !interpolation method interpolationMethodKc = IniReadInt ('interpolation', etini, section = 'crop-coefficient') IF (interpolationMethodKc == 0) THEN !data are stored in net-cdf file IF ( KeyIsPresent ('variable', etini, section = 'crop-coefficient') ) THEN kc_map % var_name = IniReadString ('variable', etini, section = 'crop-coefficient') ELSE IF (KeyIsPresent ('standard_name', etini, section = 'crop-coefficient') ) THEN kc_map % standard_name = IniReadString ('standard_name', etini, section = 'crop-coefficient') ELSE CALL Catch ('error', 'Evapotranspiration', & 'standard_name or variable missing in section crop-coefficient' ) END IF !Get the dt of imported field. Assume dt is regular dtKc = GetDtGrid (filenameKc, .TRUE.) ELSE !open file containing local measurements fileunitKc = GetUnit () OPEN ( fileunitKc, file = filenameKc (1:LEN_TRIM(filenameKc)), status = 'old') !populate metadata CALL ReadMetadata (fileunitKc, kc_stations) !set dt dtKc = kc_stations % timeIncrement !read map of crop coefficient code CALL GridByIni (ini = etini, grid = kc_code_map, & section = 'crop-coefficient', subsection = 'code-map') END IF !set time tnewKc = time END IF !set time tnewET = time ! Configuration terminated. Deallocate ini database CALL IniClose (etini) RETURN END SUBROUTINE ETinit !============================================================================== !| Description: ! Compute potential evapotranspiration SUBROUTINE PETupdate & ! ( time ) IMPLICIT NONE !Arguments with intent(in): TYPE (DateTime), INTENT(IN) :: time !local declarations: INTEGER (KIND = short) :: i, j INTEGER (KIND = short) :: jDay !!julian day REAL (KIND = float) :: sunsetHourAngle !!sunset hour angle [radians] REAL (KIND = float) :: relDist !!relative distance between the earth and the sun REAL (KIND = float) :: declination !! declination REAL (KIND = float) :: etRad !!extra terrestrial solar radiation [mm/day] !------------------------------end of declarations----------------------------- !Update Kc IF (dtKc > 0 ) THEN IF( time == tNewKc ) THEN CALL KcUpdate (time) tNewKc = tNewKc + dtKc END IF END IF IF ( time == tnewET ) THEN !it's time to update ET !set next time to update PET tnewET = tnewET + dtET IF ( compute_hargreaves ) THEN !compute values for applying Hargreaves !calculate solar declination jDay = DayOfYear(time,'noleap') declination = 0.4093 * SIN( 2. * pi * jDay / 365. - 1.405 ) !calculate sunset hour angle sunsetHourAngle = ACOS( - TAN (latCentroid) * TAN (declination) ) !calculate relative distance between the earth and the sun relDist = 1. + 0.033 * COS(2. * pi * jDay / 365.) !calculate extra terrestrial solar radiation (Allen et al. 1998) etRad = 15.392 * relDist * ( sunsetHourAngle * SIN(latCentroid) * & SIN(declination) + COS(latCentroid) * & COS(declination) * SIN(sunsetHourAngle) ) END IF DO i = 1, model_map % idim DO j = 1, model_map % jdim IF ( model_map % mat (i,j) /= model_map % nodata) THEN SELECT CASE (model_map % mat (i,j) ) CASE (HARGREAVES) CALL HargreavesSamani (time, temperatureAirDailyMean % mat(i,j), & temperatureAirDailyMax % mat(i,j), & temperatureAirDailyMin % mat (i,j), & etRad, pet % mat(i,j) ) pe % mat (i,j) = pet % mat(i,j) pt % mat (i,j) = pet % mat(i,j) CASE (HARGREAVES_MOD) CALL HargreavesSamaniModified (time, & temperatureAirDailyMean % mat(i,j), & temperatureAirDailyMax % mat(i,j), & temperatureAirDailyMin % mat (i,j), & etRad, dem % mat (i,j), & pet % mat(i,j) ) pe % mat (i,j) = pet % mat(i,j) pt % mat (i,j) = pet % mat(i,j) CASE (PRIESTLEY_TAYLOR) CALL PriestleyTaylor (temperatureAir % mat (i,j), & netRadiation % mat (i,j), & fvcover % mat (i,j), & dem % mat (i,j), & pt % mat (i,j), & pe % mat (i,j), & pet % mat (i,j) ) CASE (PENMAN_MONTEITH) !CALL PetPM (R = netRadiation % mat (i,j), & ! T = temperatureAir % mat (i,j), & ! W = windSpeed % mat (i,j), & ! rs = rsMin % mat (i,j), & ! v = fvcover % mat (i,j), & ! lai = lai % mat (i,j), & ! z = plantsHeight % mat (i,j), & ! um = relHumidityAir % mat (i,j), & ! quota = dem % mat (i,j), & ! pet = pet % mat (i,j), & ! pe = pe % mat (i,j), & ! pt = pt % mat (i,j)) CALL PenmanMonteith (temperatureAir % mat (i,j), & netRadiation % mat (i,j), & relHumidityAir % mat (i,j), & windSpeed % mat (i,j), & fvcover % mat (i,j), & plantsHeight % mat (i,j), & dem % mat (i,j), & anemometersWS % offsetZ, & hygrometers % offsetZ, & lai % mat (i,j), & rsMin % mat (i,j), & pt % mat (i,j), & pe % mat (i,j), & pet % mat (i,j) ) CASE (FAO56_PENMAN_MONTEITH) CALL FAO56PenmanMonteith (temperatureAir % mat (i,j), & netRadiation % mat (i,j), & netRadiationFAO % mat (i,j), & relHumidityAir % mat (i,j), & windSpeed % mat (i,j), & fvcover % mat (i,j), & dem % mat (i,j), & anemometersWS % offsetZ, & hygrometers % offsetZ, & lai % mat (i,j), & pt % mat (i,j), & pe % mat (i,j), & pet % mat (i,j) ) END SELECT IF ( useCropCoefficient) THEN pt % mat (i,j) = pet % mat(i,j) * kc_map % mat (i,j) !update pet pet % mat (i,j) = fvcover % mat (i,j) * pt % mat (i,j) + & ( 1. - fvcover % mat (i,j) ) * pe % mat (i,j) ELSE pt % mat (i,j) = pet % mat(i,j) END IF !check negative values IF ( pt % mat (i,j) < 0.) pt % mat (i,j) = 0. IF ( pe % mat (i,j) < 0.) pe % mat (i,j) = 0. IF ( pet % mat (i,j) < 0.) pet % mat (i,j) = 0. END IF END DO END DO END IF RETURN END SUBROUTINE PETupdate !============================================================================== !| Description: ! read new Kc values SUBROUTINE KcUpdate & ! ( time ) IMPLICIT NONE !Arguments with intent(in): TYPE (DateTime), INTENT(IN) :: time !local declarations: INTEGER (KIND = short) :: i, j, k !------------------------------end of declarations----------------------------- SELECT CASE (interpolationMethodKc) CASE (0) !update map from netcdf CALL ReadField (filenameKc, time, dtKc, dtKc, 'M', kc_map) CASE (1) !update kc map from map of crop coefficient code !read new data CALL ReadData (network = kc_stations, fileunit = fileunitKc, & time = time, aggr_time = dtKc, & aggr_type = 'ave') !interpolation DO i = 1, kc_map % idim DO j = 1, kc_map % jdim IF ( kc_map % mat(i,j) == kc_map % nodata ) CYCLE DO k = 1, kc_stations % countObs IF ( ToString ( kc_code_map % mat(i,j) ) == & kc_stations % obs(k) % id ) THEN kc_map % mat(i,j) = kc_stations % obs(k) % value EXIT END IF END DO END DO END DO END SELECT RETURN END SUBROUTINE KcUpdate !============================================================================== !| Description: ! Compute potential evapotranspiration with Hargreaves Samani model ! at daily time scale ! ! References: ! Hargreaves, G.H., and Z.A. Samani (1982). Estimating potential ! evapotranspiration. J. Irrig. Drain. Eng., ASCE, 108(3), 223-230 ! ! Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. Crop ! Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO ! Irrigation and Drainage Paper 56, 9th ed.; Food and Agriculture ! Organization of the United Nations: Rome, Italy, 1998; ISBN 92-5-104219-5. SUBROUTINE HargreavesSamani & ! (time, Tavg, Tmax, Tmin, etrad, pet) IMPLICIT NONE !Arguments with intent(in): TYPE(DateTime) , INTENT(in) :: time !!simulation time REAL (KIND = float), INTENT(in) :: Tavg !!average daily temperature [°C] REAL (KIND = float), INTENT(in) :: Tmax !!maximum daily temperature [°C] REAL (KIND = float), INTENT(in) :: Tmin !!minimum daily temperature [°C] REAL (KIND = float), INTENT(in) :: etRad !!extra terrestrial solar radiation [mm/day] !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration [m/s] !local declarations. REAL (KIND = float), PARAMETER :: HC = 0.0023 !!empirical coefficient REAL (KIND = float), PARAMETER :: HE = 0.5 !!empirical exponent REAL (KIND = float), PARAMETER :: HT = 17.8 !!empirical temperature coefficient !-----------------------------------------end of declarations------------------ !calculate evapotranspiration IF ( ( Tmax - Tmin) < 0.) THEN CALL Catch ('warning', 'Evapotranspiration', '(Tmax - Tmin) < 0.' ) pet = 0. RETURN ELSE pet = HC * etRad * (Tmax - Tmin) ** HE * (Tavg + HT) END IF !conversion from mm/day to m/s pet = pet * millimeter / day RETURN END SUBROUTINE HargreavesSamani !============================================================================== !| Description: ! Compute potential evapotranspiration at daily time scale with ! Hargreaves Samani model modified by Ravazzani et al. (2012) ! to include an elevation correction factor ! ! References: ! Hargreaves, G.H., and Z.A. Samani (1982). Estimating potential ! evapotranspiration. J. Irrig. Drain. Eng., ASCE, 108(3), 223-230 ! ! Ravazzani, G., Corbari, C., Morella, S., Gianoli, P., & ! Mancini, M.. (2012). Modified Hargreaves-Samani equation ! for the assessment of reference evapotranspiration in ! Alpine river basins. Journal of irrigation and drainage ! engineering, 138(7), 592–599. SUBROUTINE HargreavesSamaniModified & ! (time, Tavg, Tmax, Tmin, etrad, elevation, pet) IMPLICIT NONE !Arguments with intent(in): TYPE(DateTime) , INTENT(in) :: time !!simulation time REAL (KIND = float), INTENT(in) :: Tavg !!average daily temperature [°C] REAL (KIND = float), INTENT(in) :: Tmax !!maximum daily temperature [°C] REAL (KIND = float), INTENT(in) :: Tmin !!minimum daily temperature [°C] REAL (KIND = float), INTENT(in) :: etRad !!extra terrestrial solar radiation [mm/day] REAL (KIND = float), INTENT(in) :: elevation !!elevation above sea level [m] !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration [m/s] !local declarations. REAL (KIND = float), PARAMETER :: HC = 0.0023 !!empirical coefficient REAL (KIND = float), PARAMETER :: HE = 0.5 !!empirical exponent REAL (KIND = float), PARAMETER :: HT = 17.8 !!empirical temperature coefficient REAL (KIND = float), PARAMETER :: c0 = 0.817 !!empirical coefficient for elevation correlation REAL (KIND = float), PARAMETER :: c1 = 0.00022 !!empirical coefficient for elevation correlation !-----------------------------------------end of declarations------------------ !calculate evapotranspiration IF ( ( Tmax - Tmin) < 0.) THEN CALL Catch ('warning', 'Evapotranspiration', '(Tmax - Tmin) < 0.' ) pet = 0. RETURN ELSE pet = (c0 + c1 * elevation) * HC * etRad * (Tmax - Tmin) ** HE * (Tavg + HT) END IF !conversion from mm/day to m/s pet = pet * millimeter / day RETURN END SUBROUTINE HargreavesSamaniModified !============================================================================== !| Description: ! Compute potential evapotranspiration with Priestley-Taylor model ! ! References: ! Priestley C.H.B., Taylor R.J., 1972. On the assessment ! of surface heat flux and evaporation using largescale ! parameters. Monthly Weather Review, 100:81-92. ! ! Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. Crop ! Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO ! Irrigation and Drainage Paper 56, 9th ed.; Food and Agriculture ! Organization of the United Nations: Rome, Italy, 1998; ISBN 92-5-104219-5. ! ! Anderson, M. C., J. M. Norman, J. R. Mecikalski, J. P. Otkin, and W. P. ! Kustas (2007), A climatological study of evapotranspiration and moisture ! stress across the continental U.S. based on the thermal remote sensing: I. ! Model formulation, J. Geophys. Res., 112, D10117, doi:10.1029/2006JD007506 ! ! ASCE–EWRI. (2005). “The ASCE standardized reference evapotranspiration ! equation.” ASCE–EWRI Standardization of Reference Evapotranspiration ! Task Committe Rep., ASCE Reston, Va ! SUBROUTINE PriestleyTaylor & ! (airTemp, netRad, fc, elevation, pt, pe, pet) IMPLICIT NONE !Arguments with intent(in): REAL (KIND = float), INTENT(in) :: airTemp !!air temperature [°C] REAL (KIND = float), INTENT(in) :: netRad !! net radiawion [W/m2] REAL (KIND = float), INTENT(in) :: fc !! fractional coverage by vegetation [0-1] REAL (KIND = float), INTENT(in) :: elevation !! terrain elevation [m a.s.l.] !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pt !! potential transpiration (from vegetation) [m/s] REAL (KIND = float), INTENT(out) :: pe !! potential evaporation (from water or saturated soil) [m/s] REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration [m/s] !local declarations. REAL (KIND = float), PARAMETER :: alpha = 1.26 !! priestely taylor advection coefficient REAL (KIND = float), PARAMETER :: toMJm2 = 0.000001 !! factor to convert W/me to MJ m-2 s-1 REAL (KIND = float), PARAMETER :: lambda = 2.453 !2.2647 !!latent heat of vaporization (MJ / kg) REAL (KIND = float) :: netRadMJ !!net radiation flux in MJ / m2 /s REAL (KIND = float) :: es !!saturation vapor pressure REAL (KIND = float) :: des !!the slope of the relationship between saturation vapour pressure and ] REAL (KIND = float) :: airPress !!air pressure [KPa] REAL (KIND = float) :: gamma !!psychrometric constant [kPa °C-1] REAL (KIND = float) :: groundHeat !!ground heat flux [MJ m-2 s-1] !-----------------------------------------end of declarations------------------ ! compute saturation vapor pressure es = 0.6108 * EXP ( ( 17.27 * airTemp ) / ( airTemp + 237.3 ) ) !compute the slope of saturation vapour pressure curve (Allen. et al. 1998) des = ( 4098. * es ) / ( ( airTemp + 237.3 )**2. ) !compute atmospheric pressure (ASCE-EWRI, 2005) airPress = 101.3 * ( ( 293. - 0.0065 * elevation ) / 293. )** 5.26 !compute psychrometric constant gamma = 0.665 * 0.001 * airPress !compute neat radiation in MJ / m2 / s netRadMJ = netRad * toMJm2 ! compute ground heat flux as a fraction of net radiation Anderson et al. (2007) ! groundHeat is assumed to affect evaporation from bare soil groundHeat = 0.1 * netRadMJ !compute potential evaporation (from saturated bare soil or water) pe = alpha * ( des / ( des + gamma ) ) * ( netRadMJ - groundHeat) / lambda * millimeter ! (m/s) !compute potential transpiration from vegetation pt = alpha * ( des / ( des + gamma ) ) * netRadMJ / lambda * millimeter ! (m/s) !compute potential as a weighted average of bare soil and vegetated fluxes pet = ( 1. - fc) * pe + fc * pt RETURN END SUBROUTINE PriestleyTaylor !============================================================================== !| Description: ! Compute potential evapotranspiration with Penman-Monteith model ! ! Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. Crop ! Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO ! Irrigation and Drainage Paper 56, 9th ed.; Food and Agriculture ! Organization of the United Nations: Rome, Italy, 1998; ISBN 92-5-104219-5. ! ! ASCE–EWRI. 2005. “The ASCE standardized reference evapotranspiration ! equation.” ASCE–EWRI Standardization of Reference Evapotranspiration ! Task Committe Rep., ASCE Reston, Va ! ! https://wetlandscapes.github.io/blog/blog/penman-monteith-and-priestley-taylor/ ! SUBROUTINE PenmanMonteith & ! (airTemp, netRad, rh, wind, fc, z, elevation, zws, zrhum, lai, srmin, pt, pe, pet) IMPLICIT NONE !Arguments with intent(in): REAL (KIND = float), INTENT(in) :: airTemp !!air temperature [°C] REAL (KIND = float), INTENT(in) :: netRad !! net radiawion [W/m2] REAL (KIND = float), INTENT(in) :: rh !! air relative humidity [0-100] REAL (KIND = float), INTENT(in) :: wind !! wind speed [m/s] REAL (KIND = float), INTENT(in) :: fc !! fractional coverage by vegetation [0-1] REAL (KIND = float), INTENT(in) :: z !! vegetation height [m] REAL (KIND = float), INTENT(in) :: elevation !! terrain elevation [m a.s.l.] REAL (KIND = float), INTENT(in) :: zws !! wind speed measurement heigth [m] REAL (KIND = float), INTENT(in) :: zrhum !! relative humidity measurement heigth [m] REAL (KIND = float), INTENT(in) :: lai !! leaf area index [m2/m2] REAL (KIND = float), INTENT(in) :: srmin !! minimum stomatal resistance !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pt !! potential transpiration (from vegetation) [m/s] REAL (KIND = float), INTENT(out) :: pe !! potential evaporation (from water or saturated soil) [m/s] REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration [m/s] !local declarations: REAL (KIND = float) :: groundHeat !!ground heat flux [MJ m-2 s-1] REAL (KIND = float) :: es !!saturation vapor pressure (Pa) REAL (KIND = float) :: ea !!actual vapor pressure (Pa) REAL (KIND = float) :: des !!the slope of the relationship between saturation vapour pressure and temperature (kPa°C-1) REAL (KIND = float) :: airPress !!air pressure [KPa] REAL (KIND = float) :: gamma !!psychrometric constant [kPa °C-1] REAL (KIND = float) :: ws !!wind speed (m/s) REAL (KIND = float) :: z1 !!bare soil heigth (m) REAL (KIND = float) :: zm1 !! wind speed measurement heigth + z1 (m) REAL (KIND = float) :: zh1 !! relative humidity measurement heigth + z1 (m) REAL (KIND = float) :: zd1 !! bare soil zero plane displacement height (m) REAL (KIND = float) :: zrw1 !! Roughness length governing momentum transfer above zd1 (m) REAL (KIND = float) :: zrh1 !! Roughness length governing transfer of heat and vapor above zd1 (m) REAL (KIND = float) :: rabs !! aerodynamic resistance of bare soil (s/m) REAL (KIND = float) :: zm !! wind speed measurement heigth + z (m) REAL (KIND = float) :: zd !! vegetated soil zero plane displacement height (m) REAL (KIND = float) :: zrw !! Roughness length governing momentum transfer above zd (m) REAL (KIND = float) :: zh !! relative humidity measurement heigth + z (m) REAL (KIND = float) :: zrh !! Roughness length governing transfer of heat and vapor above zd (m) REAL (KIND = float) :: ra !! aerodynamic resistance of vegetated soil (s/m) REAL (KIND = float) :: cp !! humid air specific heat (MJkg^-1°C^-1) REAL (KIND = float) :: Tkv !! virtual air temperature (K) REAL (KIND = float) :: rhoAir !!air density (kgm^-3) REAL (KIND = float) :: netRadMJ !! net radiation in MJm^-2s^-1 REAL (KIND = float) :: rc !! total vegetation resistance (s/m) REAL (KIND = float), PARAMETER :: k = 0.41 !!von Karman’s constant REAL (KIND = float), PARAMETER :: epsilon = 0.622 ! ratio molecular weight of water vapour/dry air REAL (KIND = float), PARAMETER :: lambda = 2.453 !2.2647 !!latent heat of vaporization (MJ / kg) REAL (KIND = float), PARAMETER :: R1 = 0.287 !!gas specific constant (kJkg^-1K^-1) !-----------------------------------------end of declarations------------------ ! compute saturation vapor pressure es = 0.6108 * EXP ( ( 17.27 * airTemp ) / ( airTemp + 237.3 ) ) !actual vapor pressure ea = rh / 100. * es !compute the slope of saturation vapour pressure curve (Allen. et al. 1998) des = ( 4098. * es ) / ( ( airTemp + 237.3 )**2. ) !compute atmospheric pressure (ASCE-EWRI, 2005) airPress = 101.3 * ( ( 293. - 0.0065 * elevation ) / 293. )** 5.26 !compute psychrometric constant gamma = 0.665 * 0.001 * airPress !check wind speed: set minimum to 0.01 m/s IF ( wind <= 0 ) THEN ws = 0.01 ELSE ws = wind END IF !aerodynamic variables of bare soil z1 = 0.1 zm1 = zws + z1 zd1 = 2. / 3. * z1 zrw1 = 0.123 * z1 zh1 = zrhum + z1 zrh1 = 0.1 * zrw1 !compute aerodynamic resistance of bare soil rabs = LOG ( (zm1 - zd1) / zrw1 ) * LOG ( (zh1 - zd1) / zrh1 ) / ( k**2 * ws ) !aerodynamic variables of vegetated soil zm = zws + z zd = 2. / 3. * z zrw = 0.123 * z zh = zrhum + z zrh = 0.1 * zrw !compute aerodynamic resistance of bare soil IF (z == 0 ) THEN !vegetation not present ra = rabs ELSE ra = LOG ( ( zm - zd ) / zrw ) * LOG ( ( zh - zd ) / zrh ) / ( k**2 * ws) END IF ! humid air specific heat (MJkg^-1°C^-1) cp = gamma * epsilon * lambda / airPress !virtual air temperature (K) Tkv = 1.01 * ( airTemp + 273. ) !air density (kgm^-3) rhoAir = airPress / ( Tkv * R1 ) !check and convert net radiation IF ( netRad < 0 ) THEN netRadMJ = 0. ELSE netRadMJ = netRad * 0.000001 END IF ! compute ground heat flux as a fraction of net radiation Anderson et al. (2007) ! groundHeat is assumed to affect evaporation from bare soil !groundHeat = 0.31 * ( 1. - fc ) * netRadMJ groundHeat = 0.1 * netRadMJ !compute potential evaporation (from saturated bare soil or water) pe = ( des * ( netRadMJ - groundHeat ) + rhoAir * cp * ( es - ea ) / rabs ) / ( lambda * ( des + gamma ) ) * millimeter !compute potential transpiration from vegetation IF (fc > 0.) THEN IF ( lai == 0.) THEN rc = srmin / 1. ELSE rc = srmin / lai END IF pt = ( des * netRadMJ + rhoAir * cp * ( es - ea ) / ra ) / ( lambda * ( des + gamma * (1. + rc / ra) ) ) * millimeter ELSE pt = 0. END IF !compute potential as a weighted average of bare soil and vegetated fluxes pet = ( 1. - fc) * pe + fc * pt !IF (isnan (pet)) then ! write(*,*) pe, pt, des, ra, rc, lai, srmin ! pause !end if RETURN END SUBROUTINE PenmanMonteith !============================================================================== !| Description: ! Compute potential evapotranspiration with FAO56 Penman-Monteith model ! for a reference grass ! "A hypothetical reference crop with an assumed crop height of 0.12 m, ! a fixed surface resistance of 70 s m-1 and an albedo of 0.23." ! ! References: ! ! Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. Crop ! Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO ! Irrigation and Drainage Paper 56, 9th ed.; Food and Agriculture ! Organization of the United Nations: Rome, Italy, 1998; ISBN 92-5-104219-5. ! ! Lincoln Zotarelli, Michael D. Dukes, Consuelo C. Romero, ! Kati W. Migliaccio, and Kelly T. Morgan, Step by Step Calculation ! of the Penman-Monteith Evapotranspiration (FAO-56 Method), ! Agricultural and Biological Engineering Department, University of Florida, ! 2009. https://edis.ifas.ufl.edu/pdf/AE/AE45900.pdf ! ! https://wetlandscapes.github.io/blog/blog/penman-monteith-and-priestley-taylor/ ! SUBROUTINE FAO56PenmanMonteith & ! (airTemp, netRad, netRadFAO, rh, wind, fc, elevation, zws, zrhum, lai, pt, pe, pet) IMPLICIT NONE !Arguments with intent(in): REAL (KIND = float), INTENT(in) :: airTemp !!air temperature [°C] REAL (KIND = float), INTENT(in) :: netRad !! net radiation [W/m2] REAL (KIND = float), INTENT(in) :: netRadFAO !! net radiation for reference vegetation, with albedo = 0.23 [W/m2] REAL (KIND = float), INTENT(in) :: rh !! air relative humidity [0-100] REAL (KIND = float), INTENT(in) :: wind !! wind speed [m/s] REAL (KIND = float), INTENT(in) :: fc !! fractional coverage by vegetation [0-1] REAL (KIND = float), INTENT(in) :: elevation !! terrain elevation [m a.s.l.] REAL (KIND = float), INTENT(in) :: zws !! wind speed measurement heigth [m] REAL (KIND = float), INTENT(in) :: zrhum !! relative humidity measurement heigth [m] REAL (KIND = float), INTENT(in) :: lai !! leaf area index [m2/m2] !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pt !! potential transpiration (from vegetation) [m/s] REAL (KIND = float), INTENT(out) :: pe !! potential evaporation (from water or saturated soil) [m/s] REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration [m/s] !local declarations: REAL (KIND = float) :: z = 0.12 !! reference grass height [m] REAL (KIND = float) :: srmin = 70 !! minimum stomatal resistance of reference grass [s/m] REAL (KIND = float) :: groundHeat !!ground heat flux [MJ m-2 s-1] REAL (KIND = float) :: es !!saturation vapor pressure (Pa) REAL (KIND = float) :: ea !!actual vapor pressure (Pa) REAL (KIND = float) :: des !!the slope of the relationship between saturation vapour pressure and temperature REAL (KIND = float) :: airPress !!air pressure [KPa] REAL (KIND = float) :: gamma !!psychrometric constant [kPa °C-1] REAL (KIND = float) :: ws2m !!wind speed at 2m (m/s) REAL (KIND = float) :: z1 !!bare soil heigth (m) REAL (KIND = float) :: zm1 !! wind speed measurement heigth + z1 (m) REAL (KIND = float) :: zh1 !! relative humidity measurement heigth + z1 (m) REAL (KIND = float) :: zd1 !! bare soil zero plane displacement height (m) REAL (KIND = float) :: zrw1 !! Roughness length governing momentum transfer (m) REAL (KIND = float) :: zrh1 !! Roughness length governing transfer of heat and vapor above zd1 (m) REAL (KIND = float) :: rabs !! aerodynamic resistance of bare soil (s/m) REAL (KIND = float) :: cp !! humid air specific heat (MJkg^-1°C^-1) REAL (KIND = float) :: Tkv !! virtual air temperature (K) REAL (KIND = float) :: rhoAir !!air density (kgm^-3) REAL (KIND = float) :: netRadMJ !! net radiation in MJm^-2s^-1 REAL (KIND = float) :: netRadFAO_MJ !! FAO net radiation in MJm^-2s^-1 REAL (KIND = float), PARAMETER :: k = 0.41 !!von Karman’s constant REAL (KIND = float), PARAMETER :: epsilon = 0.622 ! ratio molecular weight of water vapour/dry air REAL (KIND = float), PARAMETER :: lambda = 2.453 !2.2647 !!latent heat of vaporization (MJ / kg) REAL (KIND = float), PARAMETER :: R1 = 0.287 !!gas specific constant (kJkg^-1K^-1) !-----------------------------------------end of declarations------------------ ! compute saturation vapor pressure es = 0.6108 * EXP ( ( 17.27 * airTemp ) / ( airTemp + 237.3 ) ) !actual vapor pressure ea = rh / 100. * es !compute the slope of saturation vapour pressure curve (Allen. et al. 1998) des = ( 4098. * es ) / ( ( airTemp + 237.3 )**2. ) !compute atmospheric pressure airPress = 101.3 * ( ( 293. - 0.0065 * elevation ) / 293. )** 5.26 !compute psychrometric constant gamma = 0.665 * 0.001 * airPress !move wind speed to 2 m reference height (Zotarelli et al., 2009) IF ( zws /= 2. ) THEN ws2m = wind ELSE ws2m = wind * 4.87 / LOG ( 67.8 * zws - 5.42 ) END IF !aerodynamic variables of bare soil z1 = 0.1 zm1 = zws + z1 zd1 = 2. / 3. * z1 zrw1 = 0.123 * z1 zh1 = zrhum + z1 zrh1 = 0.1 * zrw1 !compute aerodynamic resistance of bare soil rabs = LOG ( (zm1 - zd1) / zrw1 ) * LOG ( (zh1 - zd1) / zrh1 ) / ( k**2 * wind ) ! humid air specific heat (MJkg^-1°C^-1) cp = gamma * epsilon * lambda / airPress !virtual air temperature (K) Tkv = 1.01 * ( airTemp + 273. ) !air density (kgm^-3) rhoAir = airPress / ( Tkv * R1 ) !check and convert net radiation IF ( netRad < 0 ) THEN netRadMJ = 0. ELSE netRadMJ = netRad * 0.000001 END IF IF ( netRadFAO < 0 ) THEN netRadFAO_MJ = 0. ELSE netRadFAO_MJ = netRadFAO * 0.000001 END IF ! compute ground heat flux as a fraction of net radiation Anderson et al. (2007) ! groundHeat is assumed to affect evaporation from bare soil only groundHeat = 0.31 * netRadMJ !compute evaporation (from saturated bare soil or water) pe = ( des * ( netRadMJ - groundHeat ) + rhoAir * cp * ( es - ea ) / rabs ) / & ( lambda * ( des + gamma ) ) * millimeter !compute transpiration from vegetation IF (fc > 0.) THEN pt = ( 0.408 * des * ( netRadMJ - groundHeat ) + gamma * (900. / day ) / & ( airTemp + 273. ) * ws2m * (es - ea) ) / & (des + gamma * ( 1. + 0.34 * ws2m ) ) * millimeter ELSE pt = 0. END IF !compute evapotranspiration pet = ( 1. - fc) * pe + fc * pt RETURN END SUBROUTINE FAO56PenmanMonteith !============================================================================== !| Description: ! Compute actual evapotranspiration from ground ! by solving energy balance ! Original code written by Chiara Corbari and Jacopo Martinelli ! SUBROUTINE etEnergyBalanceGround & ! (airTemp, netRad, rh, wind, fc, z, elevation, zws, zrhum, lai, srmin, pt, pe, pet) IMPLICIT NONE !Arguments with intent(in): REAL (KIND = float), INTENT(in) :: airTemp !!air temperature (°C) REAL (KIND = float), INTENT(in) :: netRad !! net radiawion (W/m2) REAL (KIND = float), INTENT(in) :: rh !! air relative humidity (0-100) REAL (KIND = float), INTENT(in) :: wind !! wind speed(m/s) REAL (KIND = float), INTENT(in) :: fc !! fractional coverage by vegetation (0-1) REAL (KIND = float), INTENT(in) :: z !! vegetation height (m) REAL (KIND = float), INTENT(in) :: elevation !! terrain elevation (m a.s.l.) REAL (KIND = float), INTENT(in) :: zws !! wind speed measurement heigth (m) REAL (KIND = float), INTENT(in) :: zrhum !! relative humidity measurement heigth (m) REAL (KIND = float), INTENT(in) :: lai !! leaf area index (m2/m2) REAL (KIND = float), INTENT(in) :: srmin !! minimum stomatal resistance !Arguments with intent(out): REAL (KIND = float), INTENT(out) :: pt !! potential transpiration (from vegetation) (m/s) REAL (KIND = float), INTENT(out) :: pe !! potential evaporation (from water or saturated soil) (m/s) REAL (KIND = float), INTENT(out) :: pet !! potential evapotranspiration (m/s) !local declarations: REAL (KIND = float) :: groundHeat !!ground heat flux (MJ m-2 s-1) REAL (KIND = float) :: es !!saturation vapor pressure (Pa) REAL (KIND = float) :: ea !!actual vapor pressure (Pa) REAL (KIND = float) :: des !!the slope of the relationship between saturation vapour pressure and temperature REAL (KIND = float) :: airPress !!air pressure (KPa) REAL (KIND = float) :: gamma !!psychrometric constant (kPa °C-1) REAL (KIND = float) :: ws !!wind speed (m/s) REAL (KIND = float) :: z1 !!bare soil heigth (m) REAL (KIND = float) :: zm1 !! wind speed measurement heigth + z1 (m) REAL (KIND = float) :: zh1 !! relative humidity measurement heigth + z1 (m) REAL (KIND = float) :: zd1 !! bare soil zero plane displacement height (m) REAL (KIND = float) :: zrw1 !! Roughness length governing momentum transfer above zd1 (m) REAL (KIND = float) :: zrh1 !! Roughness length governing transfer of heat and vapor above zd1 (m) REAL (KIND = float) :: rabs !! aerodynamic resistance of bare soil (s/m) REAL (KIND = float) :: zm !! wind speed measurement heigth + z (m) REAL (KIND = float) :: zd !! vegetated soil zero plane displacement height (m) REAL (KIND = float) :: zrw !! Roughness length governing momentum transfer above zd (m) REAL (KIND = float) :: zh !! relative humidity measurement heigth + z (m) REAL (KIND = float) :: zrh !! Roughness length governing transfer of heat and vapor above zd (m) REAL (KIND = float) :: ra !! aerodynamic resistance of vegetated soil (s/m) REAL (KIND = float) :: cp !! humid air specific heat (MJkg^-1°C^-1) REAL (KIND = float) :: Tkv !! virtual air temperature (K) REAL (KIND = float) :: rhoAir !!air density (kgm^-3) REAL (KIND = float) :: netRadMJ !! net radiation in MJm^-2s^-1 REAL (KIND = float) :: rc !! total vegetation resistance (s/m) REAL (KIND = float), PARAMETER :: k = 0.41 !!von Karman’s constant REAL (KIND = float), PARAMETER :: epsilon = 0.622 ! ratio molecular weight of water vapour/dry air REAL (KIND = float), PARAMETER :: lambda = 2.453 !2.2647 !!latent heat of vaporization (MJ / kg) REAL (KIND = float), PARAMETER :: R1 = 0.287 !!gas specific constant (kJkg^-1K^-1) !-----------------------------------------end of declarations------------------ ! compute saturation vapor pressure es = 0.6108 * EXP ( ( 17.27 * airTemp ) / ( airTemp + 237.3 ) ) !actual vapor pressure ea = rh / 100. * es !compute the slope of saturation vapour pressure curve (Allen. et al. 1998) des = ( 4098. * es ) / ( ( airTemp + 237.3 )**2. ) !compute atmospheric pressure (ASCE-EWRI, 2005) airPress = 101.3 * ( ( 293. - 0.0065 * elevation ) / 293. )** 5.26 !compute psychrometric constant gamma = 0.665 * 0.001 * airPress !check wind speed: set minimum to 0.01 m/s IF ( wind <= 0 ) THEN ws = 0.01 ELSE ws = wind END IF !aerodynamic variables of bare soil z1 = 0.1 zm1 = zws + z1 zd1 = 2. / 3. * z1 zrw1 = 0.123 * z1 zh1 = zrhum + z1 zrh1 = 0.1 * zrw1 !compute aerodynamic resistance of bare soil rabs = LOG ( (zm1 - zd1) / zrw1 ) * LOG ( (zh1 - zd1) / zrh1 ) / ( k**2 * ws ) !aerodynamic variables of vegetated soil zm = zws + z zd = 2. / 3. * z zrw = 0.123 * z zh = zrhum + z zrh = 0.1 * zrw !compute aerodynamic resistance of bare soil IF (z == 0 ) THEN !vegetation not present ra = rabs ELSE ra = LOG ( ( zm - zd ) / zrw ) * LOG ( ( zh - zd ) / zrh ) / ( k**2 * ws) END IF ! humid air specific heat (MJkg^-1°C^-1) cp = gamma * epsilon * lambda / airPress !virtual air temperature (K) Tkv = 1.01 * ( airTemp + 273. ) !air density (kgm^-3) rhoAir = airPress / ( Tkv * R1 ) !check and convert net radiation IF ( netRad < 0 ) THEN netRadMJ = 0. ELSE netRadMJ = netRad * 0.000001 END IF ! compute ground heat flux as a fraction of net radiation Anderson et al. (2007) ! groundHeat is assumed to affect evaporation from bare soil !groundHeat = 0.31 * ( 1. - fc ) * netRadMJ groundHeat = 0.1 * netRadMJ !compute potential evaporation (from saturated bare soil or water) pe = ( des * ( netRadMJ - groundHeat ) + rhoAir * cp * ( es - ea ) / rabs ) / ( lambda * ( des + gamma ) ) * millimeter !compute potential transpiration from vegetation IF (fc > 0.) THEN IF ( lai == 0.) THEN rc = srmin / 1. ELSE rc = srmin / lai END IF pt = ( des * netRadMJ + rhoAir * cp * ( es - ea ) / ra ) / ( lambda * ( des + gamma * (1. + rc / ra) ) ) * millimeter ELSE pt = 0. END IF !compute potential as a weighted average of bare soil and vegetated fluxes pet = ( 1. - fc) * pe + fc * pt !IF (isnan (pet)) then ! write(*,*) pe, pt, des, ra, rc, lai, srmin ! pause !end if RETURN END SUBROUTINE etEnergyBalanceGround !============================================================================== !| Description: ! check properties for each ET model SUBROUTINE CheckSpecificProperties & ! ( model, ini ) IMPLICIT NONE !Arguments with intent(in): INTEGER (KIND = short), INTENT(in) :: model TYPE (IniList) , INTENT(in) :: ini !----------------------------end of declarations------------------------------- !check data specific for each ET model SELECT CASE (model) CASE ( PENMAN_MONTEITH ) IF ( dtTemperature <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air temperature' ) END IF IF ( dtRelHumidity <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air relative humidity' ) END IF IF ( dtWindSpeed <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires wind speed' ) END IF IF ( dtRadiation <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires net radiation' ) END IF IF ( .NOT. fvcoverLoaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires vegetation fraction coverage' ) END IF IF ( .NOT. laiLoaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires leaf area index' ) END IF IF ( .NOT. dem_loaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires digital elevation model' ) END IF IF ( .NOT. rsMinloaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires minimum stomatal resistance' ) END IF CASE ( FAO56_PENMAN_MONTEITH ) IF ( dtTemperature <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air temperature' ) END IF IF ( dtRelHumidity <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air relative humidity' ) END IF IF ( dtWindSpeed <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires wind speed' ) END IF IF ( dtRadiation <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires net radiation' ) END IF IF ( .NOT. fvcoverLoaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires vegetation fraction coverage' ) END IF IF ( .NOT. dem_loaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires digital elevation model' ) END IF CASE (PRIESTLEY_TAYLOR) IF ( dtTemperature <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air temperature' ) END IF IF ( dtRadiation <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires net radiation' ) END IF IF ( .NOT. fvcoverLoaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires vegetation fraction coverage' ) END IF IF ( .NOT. dem_loaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires digital elevation model' ) END IF !IF ( .NOT. plantsHeightLoaded ) THEN ! CALL Catch ('error', 'Evapotranspiration', & ! 'ET model requires vegetation height' ) !END IF CASE (HARGREAVES) compute_hargreaves = .TRUE. IF ( dtET /= day ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily dt' ) END IF IF ( dtTemperatureDailyMean <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily air temperature' ) END IF IF ( dtTemperatureDailyMin <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily minimum air temperature' ) END IF IF ( dtTemperatureDailyMax <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily maximum air temperature' ) END IF CASE (HARGREAVES_MOD) compute_hargreaves = .TRUE. IF ( dtET /= day ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily dt' ) END IF IF ( dtTemperatureDailyMean <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily air temperature' ) END IF IF ( dtTemperatureDailyMin <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily minimum air temperature' ) END IF IF ( dtTemperatureDailyMax <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires daily maximum air temperature' ) END IF IF ( .NOT. dem_loaded ) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires digital elevation model' ) END IF CASE (ENERGY_BALANCE) IF ( dtTemperature <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air temperature' ) END IF IF ( dtRelHumidity <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires air relative humidity' ) END IF IF ( dtWindSpeed <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires wind speed' ) END IF IF ( dtRadiation <= 0) THEN CALL Catch ('error', 'Evapotranspiration', & 'ET model requires net radiation' ) END IF END SELECT RETURN END SUBROUTINE CheckSpecificProperties END MODULE Evapotranspiration